perm filename FAITRG[PNT,HE] blob
sn#256551 filedate 1977-12-30 generic text, type C, neo UTF8
COMMENT ⊗ VALID 00013 PAGES
C REC PAGE DESCRIPTION
C00001 00001
C00002 00002 TITLE FAITRG
C00004 00003 ↑COSD: ENTRY TO COSINE DEGREES ROUTINE
C00007 00004 S2: SKIPGE A CHECK SIGN OF ORIGINAL ARG
C00009 00005 BEGIN ASIN
C00011 00006 BEGIN ACOS
C00012 00007 BEGIN SQRT
C00015 00008 BEGIN ATAN
C00016 00009 ↑ATAN: ENTRY TO ARCTANGENT ROUTINE
C00019 00010 BEGIN ATAN2
C00021 00011 SNCOS - COMPUTES THE SINE AND COSINE OF AN ANGLE BY TABLE LOOPUP
C00024 00012 START OF EXECUTABLE CODE
C00027 00013 SINE AND COSINE APPROXIMATION LOOK-UP TABLES
C00032 ENDMK
C⊗;
TITLE FAITRG
INTERN SNCOSR,SNCOSD,ASIN,ACOS,SQRT,ATAN,ATAN2
WAVE←←1 ;WAVE≠0 INDICATES COMPILE SIN,COS ROUTINES
IFN WAVE,<
INTERNAL SIN,COS,SIND,COSD
BEGIN SIN
;FLOATING POINT SINGLE PRECISION SINE AND COSINE FUNCTION
;IF THE ARGUMENT IS IN DEGREES, THE PROPER ENTRY POINTS ARE
;SIND AND COSD, WHILE IF THE ARGUMENT IS IN RADIANS, THE PROPER
;ENTRY POINTS ARE SIN AND COS.
;COSD CALLS SIND TO CALCULATE SIND(PI/2 + X)
;COS CALLS SIN TO CALCULATE SIN(PI/2 + X)
;SIND CALLS SIN AFTER A CONVERSION FROM DEGREES TO RADIANS
;THE ROUTINE CALCULATES SINES AFTER REDUCING THE ARGUMENT TO
;THE FIRST QUADRANT AND CHECKING THE OVERFLOW BITS TO DETERMINE
;THE QUADRANT OF THE ORIGINAL ARGUMENT
;000 - 1ST QUADRANT
;001 - 2ND QUADRANT, X=-(X-PI)
;010 - 3RD QUADRANT, X=-(X-PI)
;011 - 4TH QUADRANT, X=X-3*PI/2-PI/2
;THE ALGORITHM USES A MODIFIED TAYLOR SERIES TO CALCULATE
;THE SINE OF THE NORMALIZED ARGUMENT.
A← 1
B← 2
C← 3
P← 17
↑COSD: ;ENTRY TO COSINE DEGREES ROUTINE
MOVE A, -1(P) ;PICK UP THE ARGUMENT IN B
FADR A, [90.0] ;ADD 90 DEGREES
FDVR A, SCD1 ;CONVERT TO RADIANS
JRST S1 ;ENTER SINE ROUTINE
↑SIND: ;ENTRY TO SINE DEGREES ROUTINE
MOVE A, -1(P) ;PICK UP THE ARGUMENT IN B
FDVR A, SCD1 ;CONVERT FROM DEGREES TO RADIANS
JRST S1 ;ENTER SINE ROUTINE
↑COS: ;ENTRY TO COSINE RADIANS ROUTINE
MOVE A, -1(P) ;PICK UP THE ARGUMENT IN B
FADR A, PIOT ;ADD PI/2
JRST S1 ;ENTER THE SINE ROUTINE
↑SIN: ;ENTRY TO SINE RADIANS ROUTINE
MOVE A, -1(P) ;PICK UP THE ARGUMENT IN A
S1: MOVM B,A ;GET ABSF OF ARGUMENT
CAMG B, SP2 ;SINX = X IF X.LT.2↑-10
JRST SUBEXIT
FDV B, PIOT ;DIVIDE X BY PI/2
CAMG B, [1.0] ;IS X/(PI/2) .LT. 1.0?
JRST S2 ;YES, ARG IN 1ST QUADRANT ALREADY
MULI B, 400 ;NO, SEPARATE FRACTION AND EXP.
ASH C, -202(B) ;GET X MODULO 2PI
MOVEI B, 200 ;PREPARE FLOATING FRACTION
ROT C, 3 ;SAVE 3 BITS TO DETERMINE QUADRANT
LSHC B, 33 ;ARGUMENT NOW IN RANGE (-1,1)
FAD B, [0] ;NORMALIZE THE ARGUMENT
JUMPE C, S2 ;REDUCED TO FIRST QUAD IF BITS 00
TLCE C, 1000 ;SUBTRACT 1.0 FROM ARG IF BITS ARE
FSB B, [1.0] ;01 OR 11
TLCE C, 3000 ;CHECK FOR FIRST QUADRANT, 01
TLNN C, 3000 ;CHECK FOR THIRD QUADRANT, 10
MOVNS B ;01,10
S2: SKIPGE A ;CHECK SIGN OF ORIGINAL ARG
MOVNS B ;SIN(-X) = -SIN(X)
MOVEM B, C ;STORE REDUCED ARGUMENT
FMPR B, B ;CALCULATE X↑2
MOVE A, SC9 ;GET FIRST CONSTANT
FMP A, B ;MULTIPLY BY X↑2
FAD A, SC7 ;ADD IN NEXT CONSTANT
FMP A, B ;MULTIPLY BY X↑2
FAD A, SC5 ;ADD IN NEXT CONSTANT
FMP A, B ;MULTIPLY BY X↑2
FAD A, SC3 ;ADD IN NEXT CONSTANT
FMP A, B ;MULTIPLY BY X↑2
FAD A, PIOT ;ADD IN LAST CONSTANT
S2B: FMPR A, C ;MULTIPLY BY X
JRST SUBEXIT
SC3: 577265210372 ;-0.64596371106
SC5: 175506321276 ;0.07968967928
SC7: 606315546346 ;0.00467376557
SC9: 164475536722 ;0.00015148419
SP2: 170000000000 ;2**-10
SCD1: 206712273406
PIOT: 201622077325 ;PI/2
BEND ;SINE
>
BEGIN ASIN
;FLOATING POINT SINGLE PRECISION ARCSINE FUNCTION
;THE ARCSINE IS CALCULATED WITH THE FOLLOWING ALGORITHM:
; ASIN(X) = ATAN(X/SQRT(1-X↑2))
;THE RANGE OF DEFINITION FOR ASIN IS (-1.0,1.0)
;OTHER ARGUMENTS WILL GENERATE MEANINGLESS RESULTS
A← 1
B← 2
P← 17
↑ASIN: ;ENTRY TO ASIN ROUTINE
MOVN A, -1(P) ;SAVE THE NEGATIVE OF ARG
FMP A, -1(P) ;CALCULATE -(X↑2)
FAD A, [1.0] ;CALCULATE 1-(X↑2)
JUMPE A, ASIN1 ;WAS X EITHER -1.0 OR 1.0?
PUSH P, A ;NO, CALCULATE SQRT(1-X↑2)
PUSHJ P, SQRT ;ADDRESS OF ARGUMENT OF SQRT
MOVE B, -1(P) ;GET THE ARGUMENT BACK AGAIN
FDVM B, A ;CALCULATE X/SQRT(1-X↑2)
PUSH P, A ;CALCULATE ATAN(SQRT(1-X↑2))
PUSHJ P, ATAN ;ADDRESS FOR ATAN ROUTINE
JRST SUBEXIT
ASIN1: MOVE A, PIOT ;ANSWER IS EITHER PI/2 OR-PI/2
SKIPGE -1(P) ;WAS ORIGINAL ARGUMENT POSITIVE?
MOVNS A ;NO, GET -PI/2
JRST SUBEXIT
PIOT: 201622077325 ;PI/2
BEND ;ASIN
BEGIN ACOS
;FLOATING POINT SINGLE PRECISION ARC-COSINE FUNCTION
;ACOS(X) IS CALCULATED AS
; ACOS(X)=(PI)/2 - ASIN(X)
;THE RANGE OF DEFINITION FOR ACOS IS (-1.0,1.0)
;ANY OTHER ARGUMENT WILL PRODUCE AN UNDEFINED ANSWER
A← 1
P← 17
↑ACOS: ;ENTRY TO SUBROUTINE
MOVE A, -1(P) ;GET ARGUMENT IN A
PUSH P, A ;CALL ASIN ROUTINE
PUSHJ P, ASIN ;ADDRESS OF ARGUMENT
MOVNS A ;GET -(ASIN(X))
FAD A, PIOT ;CALCULATE (PI/2)-ASIN(X)
JRST SUBEXIT
PIOT: 201622077325
BEND ;ACOS
BEGIN SQRT
;FLOATING POINT SINGLE PRECISION SQUARE ROOT FUNCTION
;THE SQUARE ROOT OF THE ABSOLUTE VALUE OF THE ARGUMENT IS
;CALCULATED. THE ARGUMENT IS WRITTEN IN THE FORM
; X= F*(2**2B) WHERE 0<F<1
;SQRT(X) IS THEN CALCULATED AS (SQRT(X))*(2**B)
;SQRT(F) IS CALCULATED BY A LINEAR APPROXIMATION, THE NATURE
;OF WHICH DEPENDS ON WHETHER 1/4 < F < 1/2 OR 1/2 < F < 1,
;FOLLOWED BY TWO ITERATIONS OF NEWTON'S METHOD.
A← 1
B← 2
C← 3
D← 4
P← 17
↑SQRT: ;ENTRY TO SQUARE ROOT ROUTINE
MOVE B,-1(P) ;PICK UP THE ARGUMENT IN B
SQRT1: JUMPLE B,[SETZ A, ↔ JRST SUBEXIT]
ASHC B, -33 ;PUT EXPONENT IN B, FRACTION IN C
SUBI B, 201 ;SUBTRACT 201 FROM EXPONENT
ROT B, -1 ;CUT EXP IN HALF, SAVE ODD BIT
HRRZM B, D ;SAVE FOR FUTURE SCALING OF ANS
LSH B, -43 ;GET BIT SAVED BY PREVIOUS INST.
ASH C, -10 ;PUT FRACTION IN PROPER POSITION
FSC C, 177(B) ;PUT EXPONENT OF FRACT TO -1 OR 0
MOVEM C, A ;SAVE IT. 1/4 < F < 1
FMP C, S1(B) ;LINEAR FIRST APPROX,DEPENDS ON
FAD C, S2(B) ;WHETHER 1/4<F<1/2 OR 1/2<F<1.
MOVE B, A ;START NEWTONS METHOD WITH FRAC
FDV B, C ;CALCULATE X(0)/X(1)
FAD C, B ;X(1) + X(0)/X(1)
FSC C, -1 ;1/2(X(1) + X(0)/X(1))
FDV A, C ;X(0)/X(2)
FADR A, C ;X(2) + X(0)/X(2)
FSC A,@D ;SCALE ANSWER FOR NEWTON AND EXP
↑SUBEXIT: SUB P,[XWD 2,2]
JRST @2(P) ;EXIT
S1: 0.8125 ;CONSTANT, USED IF 1/4<FRAC<1/2
0.578125 ;CONSTANT, USED IF 1/2<FRAC<1
S2: 0.302734 ;CONSTANT, USED IF 1/4<FRAC<1/2
0.421875 ;CONSTANT, USED IF 1/2<FRAC<1
BEND ; SQRT
BEGIN ATAN
;FLOATING POINT SINGLE PRECISION ARCTANGENT FUNCTION
;ATAN(X) = X(B0+A1(Z+B1-A2(Z+B2-A3(Z+B3)**-1)**-1)**-1)
;WHERE Z=X↑2, IF 0<X<=1
;IF X>1, THEN ATAN(X) = PI/2 - ATAN(1/X)
;IF X>1, THEN RH(D) =-1, AND LH(D) = -SGN(X)
;IF X<1, THEN RH(D) = 0, AND LH(D) = SGN(X)
A← 1
B← 2
C← 3
D← 4
E← 5
P← 17
↑ATAN: ;ENTRY TO ARCTANGENT ROUTINE
MOVE A,-1(P) ;PICK UP THE ARGUMENT IN A
ATAN1: MOVM B, A ;GET ABSF OF ARGUMENT
CAMG B, A1 ;IF X<2↑-33, THEN RETURN WITH...
JRST SUBEXIT ;ATAN(X) = X
HLLO D, A ;SAVE SIGN, SET RH(D) = -1
CAML B, A2 ;IF A>2↑33, THEN RETURN WITH
JRST [MOVE A,PIOT ↔ JRST SUBEXIT]; ATAN(X) = PI/2
MOVSI C, 201400 ;FORM 1.0 IN C
CAMG B, C ;IS ABSF(X)>1.0?
TRZA D, -1 ;IF B .LE. 1.0, THEN RH(D) = 0
FDVM C, B ;B IS REPLACED BY 1.0/B
TLC D, (D) ;XOR SIGN WITH .G. 1.0 INDICATOR
MOVEM B, E ;SAVE THE ARGUMENT
FMP B, B ;GET B↑2
MOVE C, KB3 ;PICK UP A CONSTANT
FAD C, B ;ADD B↑2
MOVE A, KA3 ;ADD IN NEXT CONSTANT
FDVM A, C ;FORM -A3/(B↑2 + B3)
FAD C, B ;ADD B↑2 TO PARTIAL SUM
FAD C, KB2 ;ADD B2 TO PARTIAL SUM
MOVE A, KA2 ;PICK UP -A2
FDVM A, C ;DIVIDE PARTIAL SUM BY -A2
FAD C, B ;ADD B↑2 TO PARTIAL SUM
FAD C, KB1 ;ADD B1 TO PARTIAL SUM
MOVE A, KA1 ;PICK UP A1
FDV A, C ;DIVIDE PARTIAL SUM BY A1
FAD A, KB0 ;ADD B0
FMP A, E ;MULTIPLY BY ORIGINAL ARGUMENT
TRNE D, -1 ;CHECK .G. 1.0 INDICATOR
FSB A, PIOT ;ATAN(A) = -(ATAN(1/A)-PI/2)
SKIPGE D ;LH(D) = -SGN(B) IF B>1.0
MOVNS A ;NEGATE ANSWER
JRST SUBEXIT ;EXIT
A1: 145000000000 ;2**-33
A2: 233000000000 ;2**33
KB0: 176545543401 ;0.1746554388
KB1: 203660615617 ;6.762139240
KB2: 202650373270 ;3.316335425
KB3: 201562663021 ;1.448631538
KA1: 202732621643 ;3.709256262
KA2: 574071125540 ;-7.106760045
KA3: 600360700773 ;-0.2647686202
PIOT: 201622077325 ;PI/2
BEND ;ATAN
BEGIN ATAN2
;FLOATING POINT ARCTANGENT OF TWO ARGUMENTS
;RETURNS ARCTANGENT OF A/B
;THERE IS NO RESTRICTION ON THE ARGUMENTS
;THE ANSWER IS RETURNED IN ACCUMULATOR A.
A← 1
B← 2
P← 17
↑ATAN2: ;ENTRY POINT TO ATAN2 ROUTINE
MOVM A,-2(P)
MOVM B,-1(P)
CAML A,B
JRST FOO
MOVE A, -2(P) ;PICK UP FIRST ARGUMENT
FDVR A, -1(P) ;FORM A/B
PUSH P, A ;CALCULATE ATAN (A/B)
PUSHJ P, ATAN
SKIPL -1(P) ;IF B>0, SGN(ATAN2)=SGN(A)
JRST EXIT2 ;EXIT
JUMPGE A, ATAN2A ;IS A POSITIVE?
FADR A, PI. ;NO, SECOND QUADRANT, ADD PI
JRST EXIT2 ;EXIT
ATAN2A: FSB A, PI. ;YES,3RD QUADRANT, SUBTRACT PI
JRST EXIT2 ;EXIT
FOO: MOVN A,-1(P)
FDVR A,-2(P)
PUSH P,A
PUSHJ P,ATAN
SKIPG -2(P)
JRST XYZ
FADR A,PI2
JRST EXIT2
XYZ: FSB A,PI2
EXIT2: SUB P,[XWD 3,3]
JRST @3(P)
PI.: 202622077325
PI2: 201622077325
BEND ;ATAN2
;SNCOS - COMPUTES THE SINE AND COSINE OF AN ANGLE BY TABLE LOOPUP
BEGIN SNCOS
;THIS PROGRAM CALCULATES BOTH THE SINE AND THE COSINE OF A GIVEN ANGLE.
;THE RESULTING VALUES ARE ACCURATE TO WITH IN + OR - .0003 . RUN TIME IS
;APPROXIMATELY 260 MICROSECONDS PER CALL. THIS ROUTINE IS SOMEWHAT
;FASTER THAN THE SYSTEM SINE, COSINE ROUTINES SINCE A TABLE LOOK-UP
;TECHNIQUE IS EMPLOYED RATHER THAN THE STANDARD EXPANSION FORMULAS.
;THERE ARE TWO ENTRY POINTS INTO THIS ROUTINE: SNCOSR AND SNCOSD. SNCOSR
;IS USED FOR ANGLES GIVEN IN RADIANS AND SNCOSD FOR ANGLES IN DEGREES.
;SAMPLE CALLS ARE AS FOLLOWS:
;
; SNCOSR(THETA)
; OR
; SNCOSD(THETA)
;
;WHERE THETA IS THE ANGLE IN RADIANS OR DEGREES AS APPROPRIATE. AFTER
;COMPLETION OF SNCOS, THE SINE OF THETA IS RETURNED IN REGISTER 1 AND
;THE COSINE IN REGISTER 2. ALL NUMBERS ARE IN FLOATING POINT REPRESENTATION.
;DEFINITIONS
P←17
AC←←1
BC←←6
THETA←3
SIGN←4
SNF←5
CSF←2
SNR←7
CSR←10
TWOPI: 1.59155E-1 ;1/2*Pi
TRE60: 2.77777778E-3 ;1/360
; START OF EXECUTABLE CODE
↑SNCOSD:MOVE SIGN,TRE60 ;ENTRY POINT FOR THETA IN DEGREES
JRST .+2
↑SNCOSR:MOVE SIGN,TWOPI ;ENTRY POINT FOR THETA IN RADIANS
FMPR SIGN,-1(P) ;NORMALIZE THETA TO 0→1 = 0→360 DEGREES
MOVM THETA,SIGN ;GET THE ABSOLUTE VALUE OF THETA
FSC THETA,233-215000/1000 ;see KLFIX.FAI[SUB,SYS]
KIFIX THETA,THETA ;CONVERT NUMBER TO INTEGER 40000 = 360 DEG
LDB AC,[POINT 6,THETA,35] ;GET 6 LSB FOR FINE ESTIMATE OF THETA
MOVE SNF,SINF(AC) ;PICK UP THE FINE ESTIMATE OF SIN THETA
MOVE CSF,COSF(AC) ; " " " " " " COS "
LDB AC,[POINT 6,THETA,29] ;ISOLATE NEXT 6 BITS TO USE AS
; INDEX TO ROUGH ESTIMATE
MOVNI BC,-100(AC) ;GET INDEX STARTING FROM OTHER END OF LIST
TRCE THETA,10000 ;TEST IF THETA IN QUADRANT 2 OR 4
JRST Q2OR4 ;YES IT IS
MOVE SNR,RUF(AC) ;GET ROUGH ESTIMATE OF SIN THETA
MOVE CSR,RUF(BC) ; " " " " COS "
JRST CHKS
Q2OR4: MOVE SNR,RUF(BC) ; " " " " SIN "
MOVE CSR,RUF(AC) ; " " " " COS "
CHKS: TRNE THETA,20000 ;COMPLEMENT SIN-R IF IN QUAD 3 OR 4
MOVN SNR,SNR
TRCE THETA,30000 ;COMPLEMENT COS-R IF IN QUAD 2
TRNN THETA,30000 ;COMPLEMENT COS-R IF IN QUAD 3
MOVN CSR,CSR
MOVE AC,SNR ;COMPUTE SIN AND COS PRODUCTS
FMPR AC,CSF ;AC = SNR*CSF
FMPR SNR,SNF ;SNR = SNR*SNF
FMPR CSF,CSR ;CSF = CSR*CSF
FMPR CSR,SNF ;CSR = CSR*SNF
FADR AC,CSR ;SIN THETA = SNR*CSF+CSR*SNF
FSBR CSF,SNR ;COS THETA = CSR⊗CSF-SNR*SNF
JUMPGE SIGN,.+2 ;COMPLEMENT SIN THETA IF THETA NEG.
MOVN AC,AC
SUB P,[2(2)] ;INDICATE ARGUMENTS POPPED OFF
JRST @2(P) ;RETURN
; SINE AND COSINE APPROXIMATION LOOK-UP TABLES
RUF:
OCT 000000000000,173622052575,174621754574,175455244044,175621364657
OCT 175765311624,176454402033,176536102421,176617427017,176700560232
OCT 176761477141,177421072231,177451200614,177501153444,177530763235
OCT 177560420514,177607674250,177636757120,177665642000,177714315644
OCT 177742553514,177770564466,200407160746,200421726325,200434347317
OCT 200446640523,200460777000,200472777213,200504636312,200516331261
OCT 200527655117,200541026727,200552023632,200562641020,200573273713
OCT 200603541604,200613620034,200623504224,200633174024,200642464752
OCT 200651554614,200660441230,200667120322,200675370054,200703426274
OCT 200711251310,200716657276,200724046464,200731015364,200735542356
OCT 200742044101,200746121073,200751750124,200755350070,200760517676
OCT 200763436353,200766122762,200770354475,200772352537,200774114431
OCT 200775421546,200776471551,200777304167,200777661027,201400000000
SINF:
OCT 000000000000,165622077321,166622077321,167455457435,167622077321
OCT 167766517205,170455457326,170537667207,170622077052,170704306677
OCT 170766516500,171424353130,171455457005,171506562647,171537666473
OCT 171570772302,171622076072,171653201644,171704305374,171735411103
OCT 171766514567,172407710115,172424351723,172441013520,172455455301
OCT 172472117047,172506560600,172523222315,172537664015,172554325476
OCT 172570767143,172605430570,172622072177,172636533566,172653175136
OCT 172667636464,172704277772,172720741255,172735402517,172752043741
OCT 172766505135,173401463142,173407703615,173416124256,173424344703
OCT 173432565316,173441005716,173447226304,173455446655,173463667214
OCT 173472107536,173500330044,173506550336,173514770613,173523211053
OCT 173531431277,173537651504,173546071675,173554312047,173562532203
OCT 173570752322,173577172421,173605412502,173613632544
COSF:
OCT 201400000000,200777777766,200777777730,200777777650,200777777542
OCT 200777777410,200777777234,200777777034,200777776610,200777776340
OCT 200777776044,200777775525,200777775163,200777774573,200777774161
OCT 200777773524,200777773042,200777772332,200777771603,200777771025
OCT 200777770224,200777767376,200777766527,200777765632,200777764712
OCT 200777763747,200777762760,200777761745,200777760706,200777757623
OCT 200777756514,200777755364,200777754205,200777753005,200777751556
OCT 200777750306,200777747010,200777745467,200777744123,200777742533
OCT 200777741120,200777737462,200777735774,200777734267,200777732535
OCT 200777730756,200777727153,200777725326,200777723453,200777721556
OCT 200777717636,200777715671,200777713701,200777711665,200777707625
OCT 200777705541,200777703431,200777701276,200777677116,200777674715
OCT 200777672466,200777670213,200777665716,200777663374
BEND SNCOS
END